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I.  INTRODUCTION 


1.  Background 

Mach  reflection  phenomenon  is  of  interest  to  the  Army  because  of  its  occurrence  in  sim¬ 
ulations  of  blast  wave  interactions  typical  of  nuclear  explosions.  A  nuclear  explosion  releases 
a  spherical  blast  wave  that  interacts  with  the  ground  producing  Mach  reflection  patterns 
similar  to  the  patterns  produced  when  a  planar  shock  impinges  an  inclined  surface.  For 
planar  shock  and  wedge  interactions,  Mach  reflection  patterns  are  dominated  by  a  combi¬ 
nation  of  the  following  discontinuities:  an  incident  shock,  a  reflected  shock,  a  Mach  stem, 
slip  surfaces,  and  an  expansion  wave  from  the  wedge  leading  edge.  The  region  between  the 
Mach  stem  and  the  slip  surface  is  a  region  of  high  dynamic  pressure.  Typically  one  thinks 
of  the  crushing  eflect  of  the  incident  blast  wave  as  the  primary  defeat  mechanism  of  tactical 
equipment  experiencing  a  nuclear  explosion,  but  a  secondary  mechanism,  the  overturning  ef¬ 
fect  caused  by  high  tiynamic  pressure  regions,  must  also  be  considered  and  Army  equipment 
hardened  for  survivability. 

An  extensive  amount  of  experimental,  theoretical  and  computational  data  has  been 
published  (1-4)  for  the  reflection  of  planar  shocks  from  various  inclined  rigid  surfaces.  The 
wealth  of  qualitative  and  quantitative  data  available,  and  the  similarity  to  blast  wave  phe¬ 
nomenology,  makes  the  simulation  of  the  reflection  of  planar  shocks  from  wedge  surfaces 
a  good  choice  for  computer  code  verification  and  comparison.  Three  gas  dynamic  codes, 
BLAST2D,  SHARC,  and  STEALTH  were  run  for  several  planar  shock/wedge  geometries 
and  initial  conditions  and  their  results  compared. 

The  reflection  of  planar  shocks  from  a  sharp  compressive  corner  can  take  on  one  of 
four  psuedo-stationary  (self-similar)  patterns,  regular  reflection  (RR),  single  Mach  reflection 
(SMR),  complex  Mach  reflection  (CMR),  or  double  Mach  reflection  (DMR).  The  pattern 
that  results  depends  on  the  wedge  angle,  the  incident  shock  Mach  number,  and  the  gas 
equation  of  state.  Figure  1  presents  the  transition  boundaries  between  the  different  reflection 
patterns  for  air  with  a  perfect  gas  equation  of  state.  For  the  cases  presented  here,  Mach 
number  approximately  equal  to  2,  regular  reflection,  single  Mach  reflection,  and  complex 
Mach  reflection  are  the  dominant  patterns.  Schematics  of  these  patterns  are  also  shown  in 
Figure  1. 

Additionally,  a  shock  impinging  on  a  double  wedge  was  computationally  simulated  with 
BLAST2D  for  comparison  to  an  experimental  shadowgraph  and  SHARC  results.  When  the 
planar  shock  impinges  the  first  wedge,  a  direct  single  Mach  reflection  occurs.  However,  when 
the  single  Mach  reflection  system  impinges  the  second  wedge  an  inverse  Mach  reflection 
results.  The  interesting  feature  of  an  inverse  Mach  reflection  is  that  its  triple  point  propagates 
toward  the  wedge  surface,  and  eventually  reflects  from  the  wedge  surface  producing  a  complex 
wave  pattern  mainly  dot  .mated  by  a  regular  reflection  pattern. 

Takayamaand  Ben-Dor  (5)  and  Ben-Dor  et  al.(6)  published  experimental  shadowgraphs 
which  show  the  details  of  the  inverse  Mach  reflection  process.  Dawson  et  al.(7)  published  a 
SHARC  computational  simulation  of  the  the  double  wedge  case  in  which  the  major  features  of 


the  complex  wave  pattern  were  reproduced  in  a  density  contour  plot,  hut  some  of  the  slip  lines 
revealed  in  an  experimental  shadowgraph  were  not  computationally  captured.  Capturing 
gradients  through  contour  plots  requires  adequate  grid  resolution  for  the  problem  and  an 
adequate  choice  ol  grid  contour  levels.  In  addition,  experimental  shadowgraphs  are  based  on 
the  second  derivatives  of  density,  therefore,  computational  shadowgraphs,  similarly  based  on 
the  second  derivatives  of  density,  were  generated  using  spatially  and  temporally  second-order 
accurate  BLAST2D  results  to  see  if  all  the  slip  lines  could  be  revealed.  This  is  the  first  time 
to  this  author’s  knowledge  that  computational  shadowgraphs  for  wedge  problems  have  been 
published. 

2.  Objectives 


1  he  objectives  ol  this  report  are  to  compare  results  from  BLAST2D  to  results  obtained 
with  the  SI  E  ALT  1 1  and  SHAHC  codes  for  shocks  with  Mach  number  equal  to  2.12  impinging 
on  single  wedges.  Grid  resolution  effects  will  be  examined.  Contour  plots  and  pressure 
versus  time  histories  will  be  compared.  Also,  a  BLAST2D  computation  was  performed 
for  the  double  wedge  case  to  see  il  the  complex  inverse  Mach  reflection  pattern  produced 
experimentally  and  reproduced  partially  by  a  SHARC  density  contour  plot  could  be  fully 
captured  by  BLAST2D  results  and  a  computational  shadowgraph  contour  plot.  A  discussion 
of  the  BLAST2D  code  follows  but  for  more  information  on  the  STEALTH  and  SHARC  codes 
the  reader  should  consult  references  8  and  9. 


II.  THE  BLAST2D  CODE 


1.  Governing  Equations 


I’lie  governing  equations  for  the  blast  problems  presented  here  are  the  two-dimensional 
unsteady  Euler  equations,  written  in  integral  torin: 


d 

dt 


LQdV+L 


n  ■  FdS  —  0 


(1) 


l  he  integral  form  of  the  Euler  equations  can  be  rewritten  for  a  two-dimensional  gener¬ 
alized  cell  volume  as: 
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This  set  of  four  integral  equations  represents  the  conservation  of  mass,  momentum  in 
x  and  v  directions,  and  energy,  per  unit  volume  where  p  is  the  density,  p  is  the  pressure,  u 
and  v  are  the  velocities  in  the  x  (longitudinal)  and  v  (height)  directions  respectively,  and  e 
is  the  total  internal  energy  per  unit  volume: 

e  =  7  ~  TT  +  1/2  d(«2  +  >'1)  (4) 

(7  -  1) 

The  volume  fluxes  are  defined  as: 


l  f  =  l)r,u  -  -V'  (r>) 

V  =  y^u  +  x^v  (6) 

For  a  two-dimensional  cell,  the  integration  of  flux  over  the  surface  in  Equation  1  has 
been  replaced  in  Equation  2  by  an  integral  over  each  face  of  the  cell.  The  77-direction  is  taken 
as  the  body  normal  and  the  £ -direction  is  tangential  to  the  surface  of  the  body.  The  cell 
volume  and  walls  are  assumed  to  be  fixed  in  time.  The  metrics  x^,  j/1(  are  the  vector 

elements  of  the  cell  walls  and  V  is  the  volume  of  the  grid  cell. 

The  physical,  independent  variables  (x,y,t)  were  transformed  into  a  uniformly  spaced 
computational  grid  (£,7/,r)  by  a  general  transformation  of  the  form: 


r  =  t 


£  =  £(Ex,t/)  (7) 

i]  =  7 i(t,x,y) 


1'lte  transformations  were  chosen  so  that  the  grid  spacing  in  the  computational  space  is 
uniform  and  of  unit  length,  A£  =  1.  At;  =  1.  Thus,  the  uniform  equi-spaced  mesh  in  £  and 
7/  allows  the  use  of  unweighted  differencing  schemes.  As  a  result,  the  computational  code 
can  be  applied  to  a  variety  of  physical  geometries  and  grids. 

If  an  average  flux  is  defined  on  the  ceil  faces,  and  A£  and  At;  are  taken  as  unity,  the 
integral  form  of  the  Euler  equation,  Equation  2  can  be  rewritten  in  finite  volume  form  as: 


V.u 
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(8) 


where  the  indices  i  and  j  correspond  to  the  £  and  7/  directions  respectively  in  the  computa¬ 
tional  mesh. 

The  vectors  E  and  F  are  the  convective  numerical  fluxes  in  computational  space  (£,  7/,  r) 
consistent  with  the  physical  fluxes  E  and  F  in  (x,y,t).  The  vector  Q  consists  of  the  cell 
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averaged  dependent  variables.  The  integration  scheme  is  fully  implicit  if  m=n-(-l  and  is 
explicit  if  m=n.  The  variables  have  been  nondimensionalized  as  follows; 


x  = 


y  = 


r 

y 


u 

ii  =  — 


v 


P  = 


t  : 


P\C[ 

e 

Pici 


(9) 


'  =  —  t  =  — 

Pi  L 

where  L=l,  c=sound  speed,  and  subscript  1  represents  the  ambient  conditions  initially 
present  in  the  driven  section. 


2.  The  Computational  Algorithm 

The  fine  grid  compilations  were  performed  on  a  Cray  XMP/48  supercomputer  and 
the  coarse  grid  computations  on  an  IRIS-ID  workstation  by  discretizing  the  Euler  equations 
with  an  upwind.  Total  Variation  Diminishing  (TVD),  finite-volume,  implicit  scheme.  In 
previous  papers,  the  scheme  was  presented  in  detail  and  proved  to  be  well  suited  for  blast 
wave  calculations  (10,11).  For  a  complete  mathematical  description  of  the  algorithm,  the 
reader  is  referred  to  these  reports.  The  remainder  of  this  section  contains  a  description  of 
the  characteristics  oi  the  algorithm  in  general  terms. 

The  ULAST2D  algorithm  is  based  upon  Roe’s  approximate  Riemann  solver  (12)  cou¬ 
pled  with  upwind  flux  difference  splitting.  Other  approximate  Riemann  solvers  could  have 
been  used,  but  Roe’s  method  is  the  approach  recommended  by  Chakravarthy  when  com¬ 
putational  efficiency  is  important  (13).  The  Riemann  problem  is  given  two  initial  states 
and  the  dependent  variables  at  the  center  of  two  neighboring  control  volumes  for  a  finite 
volume  formulation,  compute  the  wave  families  (that  is  the  combination  of  shocks,  contacts, 
and  rarefaction  fans)  that  result  at  the  interface  of  the  control  volumes.  Roe’s  technique 
provides  a  direct  approximation  of  the  intermediate  states  in  contrast  to  a  solution  of  the 
exact  Riemann  problem  which  requires  an  iterative  process. 

Once  the  piecewise  constant  states  separated  by  the  wave  families  are  determined,  the 
flux  differences  across  each  wave  family  can  be  computed.  Upwinding  requires  that  the  flux 
differences  across  right  running  wave  families  (positive  eigenvalues)  be  used  in  the  derivative 
evaluations  of  fluxes  into  neighboring  fluid  cells  to  the  right  of  the  Riemann  solution  and 
that  the  flux  differences  across  left  running  wave  families  (negative  eigenvalues)  be  used  in 
the  derivative  evaluations  of  fluxes  into  neighboring  fluid  cells  to  the  left.  In  this  way  a 
method  of  characteristics-like  flavor  is  brought  into  the  numerical  algorithm.  Flux  difference 
splitting  with  TVD  was  used  to  achieve  second-order  accuracy  without  introducing  spurious 
oscillations  near  discontinuities. 

The  second-order  convective  flux  was  produced  by  adding  a  correction  term  to  the 
first-order  flux.  However,  in  order  to  avoid  spurious  oscillations,  the  correction  term  must 


4 


fulfill  the  c.iteria  for  the  algorithm  to  be  TVD  TVD  schemes  achieve  second-order  accu¬ 
racy  without  introducing  spurious  oscillations  near  discontinuities  by  employing  a  feedback 
mechanism-  'smart  numerical  dissipation’-  wherein  fluxes  are  compared  at  neighboring  con- 
tiol  volumes.  In  regions  of  little  change  no  numerical  dissipation  is  added  to  the  second 
order  correction  terms,  while  in  regions  of  large  change,  numerical  dissipation  is  added  to 
ensure  numerical  stability. 

During  this  process  new  extrema  are  not  created  by  the  numerical  dissipation.  TVD 
data  preserve  monotonicity;  a)  no  new  numerical  extrema  is  created  and  b)  the  absolute  value 
of  already  existing  numerical  extrema  must  not  increase.  TVD  schemes  have  been  rigorously 
proven  in  one  dimension,  however  the  extensions  to  two  and  three  dimensions  have  not  been 
mathematically  proven.  The  advantages  of  TVD  algorithms  over  older  schemes  are  that 
strong  gradients  and  complex  flow  fields  are  resolved  accurately  without  the  need  to  adjust 
arbitrary  smoothing  parameters.  The  disadvantages  of  upwind  differencing  with  TVD  are 
long  computing  times  caused  by  an  increase  in  the  number  of  arithmetic  operations  per 
integration  step  and  loss  of  programming  simplicity. 

The  conservative  nature  of  the  scheme  captures  shocks  and  other  discontinuities  au¬ 
tomatically.  The  finite  volume  philosophy  ensures  conservation  at  interior  and  boundary 
points.  The  scheme  is  made  implicit  by  linearizing  only  the  first-order  contribution  and  by 
employing  a  Newton  iteration  of  the  type  described  by  Rai  (14)  to  eliminate  any  approxima¬ 
tions  made.  The  implicit  version  of  the  scheme  requires  more  computations  per  integration 
step  than  the  explicit  version,  but  is  necessary  to  handle  the  stiff  nature  of  the  problems. 

The  discussion  in  this  section  assumed  the  reader  was  familiar  with  state-of-the-art 
computational  fluid  dynamic  terms  and  definitions.  Others  are  urged  to  read  the  reports 
referenced  in  this  section  for  a  full  mathematical  description. 


III.  GEOMETRY,  GRID,  AND  INITIAL  CONDITIONS  FOR  WEDGE 

COMPUTATIONS 

Computations  with  BLAST2D  were  performed  for  the  wedge  configurations  and  initial 
conditions  shown  in  Figures  2  and  3.  All  dimensions  given  are  in  centimeters.  Figure  2  shows 
the  geometries  for  the  single  wedge  cases  while  Figure  3  shows  the  geometry  for  the  double 
wedge  case,  f  igure  2(A)  was  run  with  a  189x172  grid  for  comparison  to  the  SHARC  code. 
Figures  2(B),  2(C),  and  2(D)  were  run  with  98x25,  46x25,  and  68x25  grids  respectively 
for  comparison  to  the  STEALTH  code.  The  shock  was  initialized  one  centimeter  in  front  of 
the  start  of  the  wedge  for  case  (A)  and  two  centimeters  in  front  of  the  leading  edge  of  the 
wedge  for  cases  (B),  (C),  and  (D).  Figure  3  was  run  with  a  200x200  grid  for  comparison  to 
the  SHARC  code  and  the  shock  was  initialized  at  x=.40  centimeters.  The  geometries  and 
grids  were  chosen  to  reproduce  simulations  already  completed  by  Heilig  (15),  Lottero  and 
Wortman  (16),  and  Dawson  et  al.  (7)  using  the  STEALTH  and  SHARC  codes. 


IV.  RESULTS  AND  DISCUSSION 
1.  30  Degree  Wedge  -  98x25  grid 

Figure  4  presents  BLAST2D  density  and  pressure  contours  superimposed  on  the  com¬ 
putational  grid  (black  lint's)  for  a  shock  with  Mach  number  equal  to  2.12  impinging  on  a 
30  degree  wedge.  The  contour  data  is  normalized  by  conditions  in  front  of  the  shock  and  is 
shown  at  171  microseconds  from  the  start,  of  the  calculation.  The  calculation  was  started  at 
time  zero  with  the  shock  one  centimeter  in  front  of  the  leading  edge  of  the  wedge.  For  this 
geometry  and  initial  shock,  a  single  Mach  reflection  results.  Even  at  this  coarse  resolution 
the  slip  line  is  resolved  in  the  density  contours,  but  is  too  diffuse  for  accurate  angle  mea¬ 
surements  and  comparison  to  von  Neumann  theory.  The  slip  line  is  defined  as  a  streamline 
where  the  flow  on  either  side  has  the  same  static  pressure  and  flow  direction,  however,  the 
density,  velocity  magnitude,  temperature,  and  other  related  functions  are  different. 

The  slip  line  is  more  clearly  seen  in  Figures  5  and  6  which  present  contour  plots  of 
Mach  number,  velocity  magnitude,  temperature,  speed  of  sound,  dynamic  pressure,  entropy, 
shadowgraph  function,  and  vorticity  magnitude.  Appendix  A  presents  equations  for  these 
functions.  Perhaps  the  most  interesting  of  these  plots  is  the  dynamic  pressure  plot.  Clearly 
the  region  between  the  Mach  stem  and  the  slip  line  is  the  region  of  highest  dynamic  pressure. 
The  Army  is  especially  interested  in  this  region  because  the  high  dynamic  pressure  can 
overturn  and  translate  vehicles  and  equipment,  diminishing  their  survivability  in  the  event 
of  a  nuclear  explosion. 

The  shadowgraph  function  is  the  second  derivative  of  density  in  each  direction  added 
together.  However,  beause  the  grid  is  coarse,  the  discontinuities  in  the  computational  shad¬ 
owgraph  are  extremely  smeared.  As  shown  by  the  189x172  grid  shadowgraph  contour  plot, 
higher  grid  resolution  produces  computational  shadowgraphs  that  are  comparable  to  the 
resolution  seen  in  experimental  shadowgraphs. 

In  simple  terms,  entropy  is  a  measure  of  the  disorder  of  the  system,  the  higher  the 
entropy,  the  greater  the  disorder.  The  highest  entropy  region,  shown  in  Figure  6,  is  the 
region  between  the  Mach  stem  and  the  slip  line.  The  vorticity  magnitude  plot,  Figure  6, 
shows  this  region  to  have  the  highest  vorticity  which  indicates  rotational  flow  present.  Flow 
regions  with  zero  vorticity  magnitude  are  irrotational.  A  major  difference  between  the  low 
resolution  contour  plots  just  examined  and  the  high  resolution  plots  discussed  in  the  next 
section  is  the  low  resolution  plots  indicate  the  entire  region  between  the  Mach  stem  and  the 
slip  line  to  be  rotational,  while  the  high  resolution  plots  show  only  the  slip  line  and  a  small 
area  behind  the  Mach  stem  as  rotational.  The  higher  resolution  case  discussed  next,  in  the 
absence  of  viscous  effects,  should  provide  more  accurate  results. 


6 


2.  30  Degree  Wedge  -  189x172  grid 


Figure  7  presents  BLAST2D  density  and  pressure  contours  superimposed  on  grid  lines 
for  a  fine  resolution  run.  Only  every  other  grid  line  is  plotted  for  clarity.  The  most  significant 
feature  of  the  finer  resolution  is  all  discontinuities  are  captured  more  accurately.  The  angles 
between  the  discontinuities  were  measured  and  are  in  excellent  agreement  with  Von  Neumann 
values  as  reported  by  Heilig  (15).  Figure  8  presents  Mach  number,  velocity  magnitude, 
temperature  and  speed  of  sound  contours.  Qualitatively,  the  plots  appear  to  produce  correct 
trends.  Figure  9  presents  dynamic  pressure,  entropy,  shadowgraph  function,  and  vorticity 
magnitude  contours.  The  dynamic  pressure  contour  shows  the  roll-up  of  the  slip  line  as  the 
highest  dynamic  pressure  region.  This  is  the  region,  discussed  earlier,  that  is  responsible  for 
overturning. 

3.  30  Degree  Wedge  -  Contour  Plot  Comparisons 

Figure  10  presents  a  comparison  of  BLAST2D  and  SHARC  contour  plots  for  the 
189x172  grid.  In  order  to  grid  the  wedge,  SHARC  is  restricted  to  island,  that  is  rect¬ 
angular  solid  cells,  or  shore  cells  which  are  rectangular  cells  split  along  either  diagonal  and 
half  of  which  is  treated  as  solid  and  the  other  half  as  fluid.  By  choosing  the  proper  aspect 
ratio  rectangles  and  using  shore  cells,  SHARC  could  model  a  smooth  wedge  surface  of  the 
required  slope  for  the  thirty  degree  wedge.  The  grid  for  BLAST2D  is  a  body  co:i formal  grid 
produced  using  an  elliptic  grid  generation  code. 

Both  codes  appear  to  resolve  the  discontinuities  present,  however,  the  slip  line  captured 
by  SHARC  is  not  as  crisply  defined.  Also,  the  slip  line  produced  by  SHARC  appears  to 
turn  into  the  Mach  stem  with  an  unrealistic  clustering  of  the  contour  lines.  Lottero  and 
Wort  man  are  investigating  the  cause  of  this  inaccuracy,  which  may  be  a  plotting  error  in 
the  postprocessing  routines.  Both  codes  show  some  contour  lines  that  are  numerical  noise 
in  front  of  the  incident  shock.  However,  the  noise  does  not  show  up  in  pressure-versus  time 
histories  and  is  not  considered  significant. 

Figure  1 1  presents  a  comparison  of  contours  generated  by  BLAST2D  and  the  STEALTH 
code  with  similar  resolution  in  the  y  direction,  however  the  x  direction  resolution  is  ambigu¬ 
ous  in  the  case  of  STEALTH  because  the  grid  points  move  for  this  Lagrangian  calculation. 
It  is  clear  that  the  contours  produced  by  STEALTH  are  not  as  crisply  defined  as  either 
BLAST2D  or  SHARC.  However,  it  should  be  noted  that  STEALTH  was  rezoning  in  this 
computation.  This  could  be  a  source  of  inaccuracy.  There  is  no  reason  to  believe  that 
STEALTH  numerics  alone  are  responsible  for  the  poorly  defined  results.  One  would  not 
normally  choose  a  Lagrangian  code  for  this  type  of  computation. 
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4.  30  Degree  Wedge  -  Pressure  versus  Time  Comparisons  and  Surface  Density 

Plots 

Figures  12-14  show  the  effects  of  grid  resolution  on  pressure  versus  time  histories  as 
well  as  a  comparison  of  BLAST2D,  STEALTH  and  SHARC  results.  Positions  A-F  are 
experimental  probe  locations  as  reported  by  Heilig  (15):  Position  A,  x=12.250,  y=-9.0672, 
Position  B,  x=12.250,  y=-4.3365,  Position  C,  x=  12.250,  y=-0. 78845,  Position  D,  x=18.500, 
y=-5. 7474,  Position  E,  x=lS. 500,  y=-2. 7488,  Position  F,  x=  18.500,  y=-0. 49977.  The  y  values 
are  subtracted  from  10,  the  overall  height  of  the  computational  domain,  to  get  the  positive 
y  positons.  The  leading  edge  of  the  wedge  was  defined  by  Helig  as  x=12.0,  y=0.0.  For  the 
SHARC  comparisons,  Lottero  and  Wort  man  (16)  defined  the  relative  positions  to  the  leading 
edge  to  be  the  same,  but  the  leading  edge  was  defined  as  x=0.0,  y=0.0.  As  expected,  the 
finer  grid  resolution  produces  sharper  discontinuities,  that  is,  discontinuities  with  less  rise 
time.  The  various  codes  seem  to  be  in  reasonable  agreement  for  the  same  resolution  with 
one  exception.  The  STEALTH  code  produces  a  rise  in  pressure  at  Position  B  that  occurs 
at  a  much  earlier  time,  100  microseconds,  than  BLAST2D  or  SHARC,  190  microseconds. 
Position  B  was  either  incorrectly  reported,  or  incorrectly  placed  for  the  STEALTH  run. 
BLAST2D  has  the  least  oscillatory  nature  at  discontinuities,  nevertheless  all  codes  are  in 
very  good  agreement  with  experimental  results  shown  as  solid  lines. 

Surface  density  results  from  BLAST2D  are  presented  in  Figure  15  for  the  coarse  and 
fine  resolutions  run.  The  fine  resolution  captures  discontinuities  with  smaller  rise  times, 
however  the  coarse  resolution  also  does  an  excellent  job  of  reproducing  the  correct  overall 
pressure  levels  and  trends.  Both  resolutions  are  in  excellent  agreement  with  Von  Neumann 
theory  for  pjpx  and  p3/p,. 


5.  45  Degree  Wedge  -  68x25  grid,  60  Degree  Wedge  -  46x25  grid 

The  45  and  60  degree  wedge  cases  were  run  using  BLAST2D  for  comparison  to  STEALTH 
results.  The  grids  were  chosen  to  closely  approximate  the  starting  grids  used  with  STEALTH. 
Figures  16  and  17  present  density  and  pressure  contour  comparisons  for  these  wedge  angles. 
The  45  degree  case  was  expected  to  produce  a  complex  Mach  reflection,  however,  only  a 
single  Mach  reflection  pattern  was  discernible  from  the  contour  plots.  A  kink  clearly  can  not 
be  seen  in  Figure  16.  Using  a  fine  grid  resolution,  and  a  second  order  Godunov  scheme  Glas 
ct  al.  (1),  also  produced  contour  plots  for  this  case  which  did  not  show  a  clearly  discernible 
kink,  however,  when  they  increased  the  Mach  number  of  the  shock  from  2.1  to  2.3  and  used 
the  same  fine  resolution,  a  kink  appeared.  In  general.  Figures  16  and  17  show  that  the  con¬ 
tour  plots  produced  by  STEALTH  are  more  noisy  than  those  produced  by  BLAST2D,  but 
more  disturbing  is  the  STEALTH  density  contour  plot  for  the  60  degree  wedge  case  showed 
an  attached  shock  at  the  corner  of  the  wedge  which  is  unrealistic.  BLAST2D  at  the  same 
resolution,  and  SHARC  on  a  higher  resolution  grid  (16),  showed  no  evidence  of  an  attached 
shock. 

Figures  18  and  19  present  BLAST2I)  and  STEALTH  snapshots  of  the  pressure  along  a 


constant  j  grid  line  for  the  45  and  60  degree  wedge  cases  respectively.  Also,  Von  Neumann 
predictions  are  plotted  for  reference.  Both  codes  appear  to  do  a  good  job  of  predicting  the 
overall  pressure  level  behind  the  incident  and  reflected  shocks.  Surface  density  plots  were 
computed  for  BLAST2D  results  and  are  presented  in  Figures  20  and  21  for  completeness. 


6.  25  and  60  Degree  Double- Wedge  -  199x199  grid 

Figure  22  presents  the  density  contour  plot  computed  using  BLAST2D  results  at  308 
microseconds  for  the  double  wedge  case,  initially,  the  incident  shock  reflects  over  the  first 
wedge  surface  as  a  single  Mach  reflection.  The  Mach  stem  created  by  the  first  wedge  reflects 
from  the  second  wedge  as  anot  her  single  Mach  reflection.  Later  the  triple  point  from  the  first 
wedge  reflection  process  and  the  triple  point  from  the  second  wedge  combine  to  form  one 
single  Mach  reflection  pattern  which  moves  toward  the  wedge  surface  and  is  termed  an  inverse 
Mach  reflection.  The  inverse  Mach  reflection  occurs  because  the  second  effective  wedge  angle 
is  such  that  it  cannot  support  stationary  or  direct  Mach  reflection.  For  stationary  Mach 
reflection  the  triple  point  moves  parallel  to  the  wedge  surface  and  for  direct  Mach  reflection 
the  triple  point  moves  away  from  the  wedge  surface. 

When  the  inverse  Mach  reflection  collides  with  the  wedge  surface,  a  regular  reflection 
pattern  results.  Figure  22  shows  the  density  contour  plot  produced  by  BLAST2D  for  the 
complex  pattern  that  results  after  the  collision  of  the  inverse  Mach  reflection  with  the  wedge 
surface.  When  this  figure  is  compared  to  the  experimental  shadowgraph  in  Figure  23,  one 
notices  three  slip  lines  in  the  experimental  shadowgraph  that  are  not  clearly  indicated  in  the 
density  contour  plot. 

Dawson  et  al.  published  (7)  a  density  contour  plot  for  the  same  case  using  SHARC 
which  also  did  not  capture  the  slip  lines.  Probably,  adequate  grid  resolution  and  contour 
levels  were  used,  but  gradient  plots  were  not  published.  When  a  density  gradient  magnitude 
contour  plot  was  computed  using  BLAST2D,  Figure  24,  two  of  the  slip  lines  are  resolved, 
and  furthermore  when  a  shadowgraph  function  contour  plot,  Figure  25,  was  produced,  all 
gradients  of  interest  were  resolved.  It  is  interesting  to  note  that  in  the  experimental  shadow¬ 
graph,  one  does  not  get  a  good  appreciation  of  the  differing  strengths  associated  wd.h  each 
of  the  discontinuities.  However,  the  computational  shadowgraph  shows  that  the  third  slip 
line  is  bounded  by  contour  levels  that  are  very  close  to  ambient  levels. 


V.  CONCLUSIONS 

To  conclude,  for  sharp  accurate  rise-times,  a  fine  grid  is  required,  while  for  pressure 
level  trends,  a  coarse  grid  is  adequate.  Also,  finer  resolution  showed  the  roll-up  of  the  slip 
line  for  the  30  degree  wedge  case  to  be  a  region  of  high  dynamic  pressure,  instead  of  the 
entire  region  between  the  Mach  stem  and  the  slip  line  as  indicated  by  the  coarser  resolution. 

Shadowgraph  function  contour  plots  should  be  used  to  resolve  slip  lines  that  are  too 
small  to  be  seen  in  density  and  density  gradient  contours.  The  shadowgraph  function  contour 
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plot  should  provide  the  best  comparison  to  experimental  shadowgraphs  provided  the  grid 
resolution  is  fine  enough  that  discontinuities  do  not  become  smeared. 

BLAST2D  captures  slip  lines  very  well  and  is  the  least  oscillatory  of  the  three  codes 
compared.  BLAS12D  has  the  advantages  that  grids  can  be  body  conformal,  and  that  it 
is  constructed  with  total  variation  diminishing  (TVD)  concepts  to  not  produce  numerical 
overshoots,  but  has  the  disadvantage  that  it  is  computationally  more  expensive  to  run  for  the 
same  grid  resolution.  However,  as  shown  by  the  30  degree  wedge  runs  BLAST2D  produces 
very  good  results  even  on  coarse  grids.  STEALTH  produced  an  unrealistic  attached  shock 
for  the  GO  degree  wedge  case  at  coarse  resolution  which  neither  BLAST2D  using  a  coarse 
resolution  or  SH  ARC  using  a  finer  grid  produced.  SHARC  has  the  disadvantages  of  numerical 
overshoots  unless  artificial  viscosity  is  turned  on  and  judiciously  chosen  for  the  problem, 
that  geometries  must  be  approximated  by  shore  and  island  cells  but,  is  computationally  less 
expensive  to  run  than  BLAS  T2D  for  the  same  resolution.  As  shown  by  the  30  degree  wedge 
comparisons,  the  numerical  overshoots  with  SHARC  are  significant  while  overshoots  with 
BLAS  I2D  are  not  discernible.  One  can  see  from  a  comparison  of  BLAST2D  and  SHARC 
that  the  trade-off  one  tries  to  optimize  when  running  one  computer  code  or  another  for  a 
particular  problem  is  degree  of  accuracy  versus  expense. 
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INITIAL  CONDITIONS:  Mt-  2.12  mo2-  J3023  kg/m*3 

p2-  45.7  kPa  mol-  .1064  kg/m»3 

Pi-  9.0  kPa  T2-  526.75  K 

R-  287  J/kg/K  T1-  294.65  K 

gamma-  1.4  U2-  472.74  m/s 


GEOMETRIES: 
dimensions  in  cm 

(A)  (B) 


Figure  2.  Wedge  Geometries  and  Initial  Conditions  for  Single  Wedge  Cases 


Figure  4.  Density  and  Pressure  Contours  Superimposed  on  !)S> 
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Figure  7.  Density  and  Pressure  Contours  Superimposed  on  189x172  Grid,  M;  =  2.12,  30 
Degree  Wedge 
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Figure  9.  Dynamic  Pressure,  Entropy,  Shadowgraph  Function,  and  Vorticity  Magnitude 
Jontours,  M/  =  2.12,  98x25  grid,  30  Degree  Wedge 
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Figure  li.  Density, Pressure  Contour  Plot  Comparisons  between  BLAST2D  and  STEALTH 
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Figure  12.  Pressure  versus  Time  Histories,  Positions  A-F,  BLAST2D  98x25  Grid, 
STEALTH  99x26  Grid 
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Figure  13.  Pressure  versus  Time  Histories,  Positions  A-F,  BLAST2D  139x172  Grid 
STEALTH  99x26  Grid 


23 


BLAST2D  -  189x172 


SHARC  -  189x172 


POSITION  A 
bcr<  t»o  ! 

.  /  J5  bar i e*o  i 


0.0  50.0  ICO. a  ISO. 3  230.0  210.0 


POSITION  8 

4S Z.  Sq'I.jfB1 

.  /  J  2  QCf [t«0 I 

....rVPftWV.... 

Shcrc 


3.3  w.3  loo.o  i«r  n  230.3  220.0 


I 

I 

- 1 

'  1 


3CS!TICN  C 

J  S  ^  .  C^~  [  *  »  p 

.  rkc.*WP. 

2ncrc 


0-0  20.0  tOO.O  150.0  200-0  ^0.0 


-os it; cu  o 

...?k°*WP _ 

ifvcrc 


103  ‘00-0  ISO. 3  230. c  250.0 


POSITION  e 

_ ■  *52  borii.dl 

-  •  /  •*  g  bqr  [•ooi 

. ji-PAUQ.... 

_ _ Shore _ 


o.o  so -o  iao.0  iso. a  200.0  250.o 

TIttC  1 •icro««c I 


POSITION  r 
.  «S7  bor <t»pl 

•  Bgn4«IaZ 

E3®E 


0.0  50.0  100. a  ISO.Q  200.0  2SO.O 

TIttC  (■Lcroitc) 


Figure  14.  Pressure  versus  Time  Histories,  Positions  A-F,  BLAST2D  189x172  Grid, 
SHARC  189x172  Grid 
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Figure  15.  Surface  Density  Plots  -  BLAST2D  with  98x25  Grid  and  189x172  Grid 
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Figure  19.  Snapshot  of  Column  35  Pressure  Comparison,  BLAST2D  and  STEALTH, 
=2.12,  60  Degree  Wedge 


Figure  21.  Surface  Density  Plot,  BLAST2D,  60  Degree  Wedg( 
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Figure  24.  Density  Gradient  Magnitude  Contour  Plot, 
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Figure  25.  Shadowgraph  Function  Contour  Plot,  25  and  60  Degree  Double  Wedge 
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APPENDIX  A:  Function  Definitions 
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p  -  static  pressure 
p  -  density 

u  -  x  component  velocity 
v  -  y  component  velocity 
Cv  -  specific  heat  at  constant  volume 
Cp  -  specific  heat  at  constant  pressure 
Velocity  magnitude 


V  =  \/u2  +  v2 

(A-l) 

Speed  of  sound 

a 

II 

(A— 2) 

Mach  number 

a 

(A— 3) 

Dynamic  pressure 

Q  =  1/2  p(u2  +  u2) 

(A-4) 

Entropy 

S  =  Cvln—  +  CJn  — 

Poo  P 

(A-5) 

Shadowgraph  function 

d2p  d2p 

1  dx 2  dy2 

(A-6) 

Vorticity  magnitude 

dv  du 

Q  =  _  _  _ 
dx  dy 

(A— 7) 
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